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We present a general theory for predicting the interaction potentials between DNA-coated colloids, and 
more broadly, any particles that interact via valence-limited ligand-receptor binding. Our theory correctly 
incorporates the configurational and combinatorial entropic factors that play a key role in valence-limited 
interactions. By rigorously enforcing self-consistency, it achieves near-quantitative accuracy with respect to 
detailed Monte Carlo calculations. With suitable approximations and in particular geometries, our theory 
reduces to previous successful treatments, which are now united in a common and extensible framework. We 
expect our tools to be useful to other researchers investigating ligand-mediated interactions. A complete and 
well-documented Python implementation is freely available at http://github.com/patvarilly/DNACC 



I. INTRODUCTION 

The theoretical study of the statistical behavior of soft 
matter systems generally begins with a description of the 
interaction between their constituent particles. One ex- 
citing system that has emerged in recent years is that of 
DNA-coated colloids (DNACCs), where interactions be- 
tween colloids are driven by partial hybridization of com- 
plementary ssDNA strands grafted on their surface^} 
To date, DNACCs have been used mainly for building 
highly sensitive biomedical probes 5 -' and to assemble bi- 
nary nanoparticle crystals of numerous morphologies^^!. 

The aim of this paper is to present a general theory of 
DNACC interactions that can be applied to any DNACC 
system, using only mild assumptions that are valid for 
common DNACCs. A key feature of our approach is that 
once the binding statistics of an arbitrary pair of DNA 
strands is known, all combinatorial and competitive ef- 
fects can be computed separately. Disentangling these 
two aspects is what gives our theory its generality. Our 
efforts build on past theoretical attempts to understand 
DNACC experiment s in t erms of the statistical mechan- 
ics of bond format iorP^HIl w hich have all partially suc- 
ceeded in modeling various experimental setups, differing 
for instance in coating densities, nucleotide sequences, 
strand lengths, stiffness and colloid sizes. To date, how- 
ever, a rigorous an d gen eral solution to the problem has 
not been presentecP^Hlll. 

A salient feature of DNACCs is that their interaction 
is mediated by surface-grafted ligands, each of which can 
form a single bond at a time. Owing to this restriction, 
we say that such interactions are "valence-limited"^. 
While we phrase our theory in the language of DNACCs, 
the treatment is equally applicable to any valence-limited 
interaction, including multi-valent particles for drug de- 
livery, viruses and cells^H 2 ^. 

The rest of the paper is structured as follows. In Sec- 
tion [nj we present a statistical mechanical derivation of 
our theory in its most general form, followed by a spe- 
cialization to DNACCs and a computationally simpler 
mean-field approximation. In Sect ion [TTH we compare the 



results of our theory with those from more accurate and 
computationally expensive Monte Carlo simulations. We 
have chosen as examples specific DNAC C systems that 
have been presented previously^ * 17 * 25 -!, as well as two 
novel ones, to show that they can all be described using a 
single, unifying framework. In Section [TV] we summarize 
our approach and reproduce the essential equations for 
reference purposes. Finally, in Section[Vj we discuss pos- 
sible applications of our theory, both for DNACCs and 
for more general valence-limited interactions. 



In Appendix |Aj we discuss an important technical as- 
pect of our theory, namely, how to calculate the entropic 
cost of binding for two general polymeric ligands. We 
also summarize the results for short, stiff rods (which 
model s hort d sDNA tethers common in micron-sized 
DNACCfpmSl) and ideal chains ( as uged in Refs. Eland 

n~5|) . Appendix [B] explains how our theory is connected 
to the Local Chemical Equilibrium treatment recently 
introduced by Rogers and Crocker^ (itself an improve- 
ment of the work of Biancaniello, Kim and Crocker in 
Ref. [3]) . We show how a physically motivated and sim- 
ple change to the equations used by these authors cor- 
rects their theory's principal shortcoming, and that the 
corrected version reduces to a specialization of our the- 
ory. Finally, Appendix [C] presents an alternate deriva- 
tion of the mean-field approximation to our theory using 
a saddle-point analysis, making explicit the connection 
between our theor y and earlier successful but more spe- 
cialized approache d 16 * 25 !. 



Although this paper explains our approach in full 
detail, we emphasize that the final results are simple 
and can be used widely by other researchers to calcu- 
late pair potentials and binding configurations of gen- 
eral DNACCs. To this end, we have made freely avail- 
able a complete and well-documented Python imple- 
mentation of the theory, which can be found at http : 
//github . com/patvarilly/DNACC| 



II. THEORETICAL MODEL 

Under experimentally relevant conditions, DNA- 
coated colloids interact through entropi c repulsion and 
hybridization-mediated attractiorPHSHS Whereas the 
repulsion is simple to estimate for common geometries, 
as shown in Appendix |Aj calculating the attraction is 
more complicated. The key subtlety is that the attrac- 
tion is mediated by valence- limited binding, i.e., each 
DNA strand can bind at most one partner at a time^. 
In this section, we first present a general treatment of 
this attraction for generic valence-limited ligands. Then, 
we specialize our treatment to DNACCs. Finally, we 
construct a mean-field approximation that is suitable for 
geometries with translational symmetry. 

A. Statistical Mechanics of valence-limited interactions 
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Ligand / Strand 



Consider a system of N generalized ligands, depicted 
in Figure [lj a). At any given time, an individual ligand i 
can bind at most one other ligand j , with a resulting free 
energy change AG^-. If i and j cannot bind, either for 
geometric or chemical reasons, we set AG^- to infinity. 
Since i cannot bind to itself, we also set AG^ to infinity 
for all i. The set of all AG^- values completely charac- 
terizes the relevant chemical detail of the system. We ex- 
plicitly assume that AG^ is independent of the presence 
of other ligands, which is accurate when the density of 
ligands is low enough that two ligands not bound to each 
other can be thought of as not interacting at all. This 
condition holds for the grafting densities and salt con- 
centrations typicall y used in DNACC experiments with 
micron-sized colloids^ 1 12 1 14 1 27 1 . 

Let 4> be a binding configuration of the system, i.e. a 
list of bound ligand pairs (i, j). Taking the state where 
no bonds are formed as a reference state, the free energy 
of (j) is the sum of the binding free energies of each bond. 
Summing over all (f>, we obtain the exact partition func- 
tion of the system, Z, relative to that of the reference 
state, Zq: 



Z 



e n «- 



0AGi 



(1) 



For a particular ligand i, the terms in the sum above 
can be classified according to the binding partner of i, 
if any. Let Z_i be the partition function for a system 



where i is unbound, and let Z_ 



h-3 



likewise be the par- 



tition function for a system where both i and j are un- 
bound. The sum of the terms in Eq. ([I]) where i is un- 
bound is Z-i/Zo, whereas the terms where i is bound 
to a ligand j out of the remaining N — 1 ligands is 
exp(— f3AGij) Z_i_j/Z$. With these definitions, it fol- 
lows from Eq. ([TJ) that 



Z-t 
Z 



-PAGi 



for all i, (2) 



FIG. 1. (a) Schematic of a DNACC system. Any pair of lig- 
ands, i and j, can form a bond at free energy cost AGij , which 
is approximately independent of the presence of the other lig- 
ands. (b) Anatomy of a typical DNACC ligand, consisting of 
a long tether capped by a point-like sticky end. Throughout, 
we use the terms ligand and strand interchangeably. 



where the sum is over the N — 1 ligands j distinct from 
i. Let pi be the probability that ligand i is unbound: 



Pi 



Z-i/Z. 



(3) 



In Eq. ([2]), the factor is the probability that 

both i and j are unbound. Noting that Z_i^_j /Z_j is the 
conditional probability for i to be unbound given that j 
is unbound, we can rewrite Z-i t _j/Z as 



Z—i—j _ Z—j-j Z_j 
Z Z-j z 



Z-i 



-z~z =1M >*' 



(4) 



Above, we have approximated the binding of i and j as 
uncorrelated events. This approximation is accurate if 
both i and j have many possible binding partners with 
comparable binding free energies. 

Using Eqs. Q-Q, we obtain the following set of N 
self-consistent equations for the quantities {pi}: 



l=Pi+J2 e ~ f}AGii PiPi> 



for all i. 



(5) 



When rewritten in the following convenient form, these 
equations can be solved, for example, by self-consistent 
iteration: 



1 



Pt = 



for all 



(6) 



The probability Pij that ligands i and j are bound is then 
given by 



Pij = e pAGi3 PiPj. 



(7) 
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Eqs. Q and Q are a principal result of this paper. It 
can be shown that these results are a single-particle ana- 
log of the more familiar chemical equilibrium conditions 
between the concentrations of N species of particles, with 
binding constants proportional to exp(— /3AGjj). 

Eq. © can be recast as a minimization problem, which 
leads to a more efficient method of solution than fixed- 
point iteration. A similar observation has recently been 
made when solving the equations that arise in umbrella 
sampling free energy calcula tions ( WH AM/MB AR), 
which have a similar structur e! 28 * 29 !. Additionally, this 
form allows us to prove that Eq. ^ has exactly one so- 
lution, which satisfies < pi < 1 for all i by construction. 
We define f({pi}), the function to be minimized, as fol- 
lows: 



f({Pi}) = E^ ~ lnpi ^ + o E ftfci ^'' 



(8) 



where kij = exp(— /3AGij). The stationary points of / 
are given by the solutions to the following system of equa- 
tions, 

^- = 1-1 

dpi p, 



= 



E 

3 



ijPj , 



(9) 



which are equivalent to Eq. ([6]). We now demonstrate 
that all stationary points of / are local minima, imply- 
ing there is a single stationary point that is the global 
minimum. To this end, we prove that the Hessian of / at 
a stationary point, d 2 f jdpi dpj, is positive definite, im- 
plying that the stationary point is a local minimum, by 
showing that the quantity J^ij Vi(d 2 f / dpi dpj)vj is pos- 
itive for all nonzero A-component vectors {vi}. Indeed, 



E* 

ij 



d 2 f 

1 dpi dpj 



E^ 

ij 
ij 



Pi 



Pi 



[5-ij + PikijPj] 



Uj, 



Pj 



(10) 



(11) 



The non-negative quantities in square brackets are the 
elements of a strictly diagonally dominant, symmetric 
matrix with positive diagonal elements, which is thus 
positive definite^!. Indeed, the i th element of its diag- 
onal is equal to 1 + pfka > 1 (we allow nonzero ku here, 
in anticipation of the mean-field self-consistent theory of 



the 



Eq. (21) below), while the sum of the off-diagonal terms 



row satisfies 



y~]pjkjjPj -Pihi < yiPihjPj = 1 ~Pi < 1) (12) 

3 3 

where the equality and the strict inequality both follow 
from Eq. ^ at the stationary point. Hence, the right- 
hand side of Eq. ( 11 ) is positive for every nontrivial choice 



of {«,}, so the Hessian of / is positive definite, as claimed. 

Having calculated the probabilities of each bond form- 
ing (Eq. ([7])), we now calculate the resulting attraction 
free energy, F att , given by 



f3F a - 



Hz/z ). 



(13) 



This quantity is most easily estimated using thermody- 
namic integration. One possible integration path is to 
offset all interaction energies by an amount A//3, i.e. 



/3AGij ^ (3AGij + A. 



(14) 



Let Z(X) be the partition function after this replacement. 
When A is infinite, Z(X) = Zq, so 



dX d[-Hz(\)/z )} 



dA ^Py(A), 

i<3 



(15) 

where Pij(X) is given by Eqs. (|6| and |7]) after effecting 
the substitution of Eq. (14). 



Once the free energy is known at one temperature, a 
more physical integration path can be used to obtain it 
at any other temperature. To proceed, we denote the 
dependence of bond energies on (3 explicitly with the no- 
tation AGij((3), and similarly for Pij((3) and F att ((3). We 
obtain 



PiFmWi) - PoF^tiPo) 



i<3 



(16) 

In calculating the integral in Eq. (fl6| from f3 to j3\, the 



free energy of binding at all intermediate temperatures 
is obtained as well. 

Once i^tt is known, the free energy of entropic repul- 
sion, -Fi-cp, must be added to obtain a full interaction free 
energy (see Appendix [A] for details). Additionally, it is 
usually convenient to offset all free energies to yield zero 
total energy when the colloids are an infinite distance 
apart, which we do in Section [III| 

In the limit of weak binding, when all the factors 
exp(— f3AGij) are very small, the above results sim- 
plify significantly, as all the complications of binding 
competition disappear. Discarding all terms of order 
exp(— 2/3AGij) and higher, we obtain 



Pij fa cxp(-f3AGij), 

/3/att ~ - E ft J- 
i<3 



(17) 

(18) 
(19) 



In particular, the last relation can be used to sidestep 



the thermodynamic integration in Eq. ( 15 ), and is related 



to a similar "P oisson" approximation used in the litera- 
tur d 3 * 10 * 13 * 17 * 19 *. However, we have previously pointed out 
that under realistic experimental conditions, this weak- 
binding approximation is rather pooi^, so its utility for 
experiments is limited. 



B. Specializing to DNA-coated colloids 

For a large class of ligands, binding occurs only in a 
small portion of the ligand. In DNACCs in particular, 
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the DNA strands that are used consist of a long, inert 
tether capped by a reactive, short sticky end, shown in 
Figure [TJb). The thermodynamics of association of DNA 
oligomers is well understood^, so for a pair of strands, 
i and j, we can predict from their sequence the value of 
the hybridization free energy of their sticky ends in so- 
lution, AG?-. Separating out this contribution to AG 



we have 



AG,; 



AG?, + AG, 



(cnf) 



(20) 



where AGy is the free energy penalty that arises from 
the reduced configurational space avail able t o two bound 
ligands with respect to two free ligands^EH. This quan- 
tity depends only on the polymer statistics of the tethers 
and the way in which they are attached to the colloids 
(e.g., tethered on the surface vs. mobile), but not on the 
sticky ends, and in many cases is a completely entropic 
effect. In Appendix |XJ we list explicit expressions for 
AG-™^ applicable to short, rigid dsDNA tethers grafted 
on plates. We also show how to use common polymer 

statistics techniques to estimate AG[" 1 ^ to arbitrary ac- 
curacy for more complicated tet hers, such as long, flexible 
ssDNAs modeled as ideal chaing2B2]. 



C. Mean field theory for DNA-coated plate interactions 

In some regimes, notably micron-sized colloids with 
short DNA tethers, the range of the binding interactions 
is small compared to the radii of curvature of all the 
colloidal surfaces involved. In such cases, we can use 
the Derjaguin approximation^ to calculate particle in- 
teraction energies, knowing only the interaction energy 
density of a pair of uniformly coated parallel plates as a 
function of plate separation. For this translationally in- 
variant geometry, a substantial simplification is possible 
if the ligands are treated at a mean field level rather than 
individually. 

The simplified system consists of two parallel plates 
separated by a distance h. We consider M different types 
of ligands, labeled by letters a and b below. Two ligand 
types can differ in their sticky end sequences, tether types 
or in the plate on which they are grafted. Correspond- 
ing to these M types, there are M(M + l)/2 potentially 
different solution hybridization free energies, denoted by 
AG° b . Let a a be the grafting density of a-type tethers, 
which has units of tethers per unit area. Since we approx- 
imated the system as uniform, all a-type tethers have the 
same probability of being unbound, denoted by p a . By 
approximating the sum in Eq. ^ by a continuous inte- 
gral, we obtain a mean-field, self-consistent set of equa- 
tions 

1 



Pa = 



1 + J2b K ab&bPb ' 



where 



J d 2 r b exp[-/?AGi: 



(cnf), 
b \ 



(21) 



(22) 



In the last equation, AG^ nf ^ is the configurational cost 
for a fixed a-type ligand to bind a 6-type tether grafted 
at rv The factor K a b has units of area, and can be 
interpreted as an average Boltzmann factor for tether 
binding times the area on plate B within reach of a fixed 
a-type tether. A relation analogous to Eq. |7]) yields the 
density of bonds between a-type and 6-type tethers, a ab : 



Vab = VaPa.K a bPb<Jb- 



(23) 



Similarly, using a relation analogous to Eq. (15), we can 



use thermodynamic integration to estimate the free en- 
ergy of attraction per unit area, f a tt(h), at fixed plate 
separation: 



f c 

tW = - / 

Jo 



dA 2 ( cr ab (\), 

a<b 



(24) 



where A is the integration parameter defined in Eq. (14 1. 

Eqs. (pT|)-(24) are a second principal result of this 
paper. They constitute a self-consistent mean field the- 
ory for the DNA-mcdiated attraction of uniformly coated 
parallel plates. In Appendix[Cj we show that, when spe- 
cialized to only one or two kinds of bonds, these equations 
reduce to the earlier self-consistent mean field theories of 
Refs. \W\ and [25] In these particular cases, we also show 



how Eqs. ( 21 1— ( 24 ) may be derived from a saddle-point 



approximation of a mean-field partition function. 

For the case of short, rigid dsDNA tethers, the fac- 
tors K a b have particularly simple expressions, derived in 
Appendix |Bj In particular, let tethers of type a and b 
have lengths L a and Lb- Without loss of generality, let 
L a < Lb- When a and b are on the same plate, 



K ab = 



Po 




L b < h; 
otherwise. 



(25a) 



whereas when a and b are on opposite plates, 



K ab = 



-/3AG" 



Po 



0. 



L a + L b < h; 



La l a L i; h , L b <h<L a + L b ; 
l/h, 



otherwise. 



(25b) 



When the tethers are modeled as freely-jointed chains 
of segments of length £, a choice that has been suggested 
for modeling ssDNA tethers^, the factors K a b can be 
approximated with the following expressions, derived in 
Appendix |A 2| Let iV a and N b be the number of segments 
in tethers of type a and 6. When a and 6 are on the same 
plate, 



K, 



a b 



-f3AG° b 



X I/ ir(N a + N b )£ 2 61 



3h 2 

+ N b )£ 2 



/-(S0-(i 



/ 3fe 2 
2N h f 



(26a) 
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whereas when a and b are on opposite plates, 

„-^ AG - r ih 2 

K a b 




(26b) 



In both of these equations, the error function, erf (a;), is 
defined as 



erf(x) = — / die 



(27) 



As before, when all the factors -Kabfb and a a Kab are 
very small, one can implement a weak-binding approx- 
imation to sidestep the thermodynamic integration in 
Eq. (15), which yields 



b 

Vab ~ GaKabO'b, 

/3/att ~ — ^ &ab- 
a<b 



(28) 

(29) 
(30) 



However, as pointed out previously, the utility of this 
approximation under realistic experimental conditions is 
limited. 

Frequently, DNACCs are covered uniformly at a fixed 
density a with a single type of strand, and only the sticky 
ends vary from one colloid to the next. In this case, we 
can obtain an explicit mean-field expression for /3/ a tt be- 
tween a plate grafted with a-type strands and one grafted 
with fe-type strands. Defining J as the dimensionless fac- 
tor 2K a b(T, we obtain 



Cab = 



J 



(i + J) - VT+2J 



(31) 



-(1 + J) + V1 + 2J 

+ 2Jln{i(l + V1 + 2J)} 



(32) 



This expression makes explicit the breakdown of the 
weak-binding approximation (in this case, /3/ a tt ~ — Gab) 
when J > 1. Unfortunately, even for the simple case 
of rod-like DNA strands, the above expression does not 
yield a simple analytic result for the pair potential be- 
tween to DNA-coated spheres via the Derjaguin approxi- 
mation. Since this difficulty only increases when we allow 
for uneven coverages and multiple potential interactions, 
we believe that our numerical thermodynamic integration 
procedure is justified. 



III. RESULTS AND DISCUSSION 

In this section, we compare the results of detailed 
Monte Carlo simulations (MC) with those derived from 
our theory, described above. We refer to the results of the 
theory where each strand is modeled individually (Eqs. 
je]), Q and |l5| )) as the SCT results, and those of the 
form where the strands are modeled at a mean field level 
(Eqs. plll-pl) as the SCTMF resultsP. 



A. Colloids coated with short dsDNA strands 

The simplest possible DNACC is a flat plate that is 
uniformly covered by a single type of DNA strand. In- 
teraction energies obtained from plate geometries can be 
used to approximate interaction energies in more com- 
plex geometries by use of the Derjaguin approximation. 
Although, our theory can be used in arbitrary geometries, 
in this section we focus on plates for clarity. 

The interaction free energy per unit area of two such 
DNACCs has be en ca lculated previously from Monte 
Carlo simulation d 15 * 16 l In this case, the tethers were 
modeled as short, stiff rods of length L = 20 nm, grafted 
at a density of l/(20nm) 2 . Figure [2] shows the fraction 
of bonds formed between the two plates when separated 
by a distance h = L, as a function of the sticky end so- 
lution hybridization free energy, AG . The results are 
calculated both from Monte Carlo simulations and from 
our theory, using an identical realization of the tether 
grafting positions on a periodic 200 nm x 200 nm plate. 
The figure inset shows the free energy per unit area of 
these plates, as a function of plate separation, for three 
different values of AG . 

For the grafting density used in Figure [2j each tether 
typically has around 9 possible binding partners, so the 
mean- field approximation is excellent. Small discrepan- 
cies arise only at low temperature and large plate sepa- 
ration, where the number of possible binding partners is 
low and Poisson fluctuations in local coverage density are 
relevant. Such fluctuations are naturally accounted for in 
the explicit self-consistent theory, whose predictions are 
thus nearly indistinguishable from MC simulations in all 
cases. 

Figure [3] shows the analogous results for a pair of 
DNACCs coated with two kinds of tethers which has 
been described previously^. We label the tethers on 
one plate by a and b, and those on the other plate 
by a' and b'. Each tether type is grafted at a density 
of 1/(20 nm) 2 . The sticky ends of a and a' hybridize 
strongly, with free energy AG , forming bridges between 
the two plates. Conversely, the sticky ends of a and 6, 
as well as those of a' and hybridize weakly, with free 
energy AG + 5 &bT, forming loops within each plate. A 
competition thus arises between strong bridges and weak 
loops, which has been shown to result in interaction po- 
tentials with an unusual temperature dependence that 
should significantly enhance the ability of these DNACCs 
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Solution hybridisation free energy AG (k^T) Solution hybridisation energy AG 



FIG. 2. Density of bonds formed between two fiat plates 
coated uniformly with short, rigid dsDNA tethers of length 
L = 20 nm, at plate separation h — L, as a function of the 
sticky end solution hybridization free energy, AG . Inset: 
potential of mean force per unit area as a function of plate 
separation, for three different values of AG . 




Solution hybridisation free energy AG (fcaT) 



FIG. 3. Density of strong and weak bonds formed between 
plates with competing-linkages, at plate separation 20 nm. 
Strong a-a' inter-plate bonds (bridges) have sticky end hy- 
bridization free energy AG , and weak a-b and a'-b' intra- 
plate bonds (loops) have sticky end hybridization free en- 
ergy AG + 5 ksT. See text for details. Inset: potential 
of mean force per unit area as a function of plate separation, 
for AG = -5k B T (upper), -15fc B T (bottom) and -25fc B T 
(middle). When the sticky end binding is strongest, the plates 
barely interact. 

to form binary crystals^. All the features of these inter- 
actions are reproduced with quantitative accuracy by our 
self-consistent theory. 

The key advantage of the present theory is its gener- 
ality and comparative simplicity. With little additional 
effort, we can consider systems with many different kinds 
of grafted DNAs, each with different entropic costs for 
binding, and with different and potentially nonuniform 



FIG. 4. Density of strong, weak and c-c' bonds formed in the 
third system, at plate separation 20 nm (see text for details) . 
Inset: potential of mean force per unit area as a function of 
plate separation, for three values of AG . 

coating densities. In the previous two cases, the result- 
ing SCTMF theory is equivalent to earlier self-consistent 
mean-field theories that were less general and that were 
derived more heuristically (Appendix [C]). In our third 
and fourth examples, we go beyond what has been done 
previously to illustrate the generality of our approach. 

For our third example, we generalize the setup of Fig- 
ure [3] to include a third type of tether on each plate, 
respectively c and c'. The sticky ends of c and d hy- 
bridize only with each other and with the same strength 
as the a and a' sticky ends. Furthermore, we set the 
lengths of the a, b, a' and b' tethers to 30 nm and those 
of the c and d tethers to 10 nm. Finally, we set the graft- 
ing density of the c and d tethers to 1/(30 nm) 2 . Thus, 
our third system features three different types of rods, 
of different lengths, with different coverages and with a 
nontrivial interaction matrix. Figure [4] summarizes the 
behavior of this system. Our choices result in the inter- 
action energy per unit area having up to two different 
minima, depending on the specific temperature (binding 
strength) . The example also illustrates one way in which 
the mean-field approach fails, and how the full theory be- 
haves for such cases. The large discrepancies in the num- 
ber of c-c' bonds formed (black solid and dashed curves 
in Figure[4|, and the resulting discrepancies in the poten- 
tials of mean force, result from the low coverage density 
of these tethers. The typical separation between graft- 
ing points that result from Poisson statistics is far larger 
than the length of the tethers, but these inhomogencitics 
are not captured by a mean field approach. Importantly, 
since the full theory retains detailed information on the 
grafting points of c and d tethers, the MC results are 
again recovered to quantitative accuracy. 

Our fourth and final example in this section is a sys- 
tem of two parallel plates, each with a 100 nm-radius 
circular patch of dsDNA strands, distributed uniformly 
within the patch at a coverage density 1/(20 nm) 2 . As 
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FIG. 5. Potential of mean force for two parallel plates coated 
with circular patches (radius 100 nm) of short dsDNA strands 
with complementary sticky ends. The coating density is 
1/(20 nm) 2 . The sticky ends have solution hybridization free 
energy of -6k B T (top), -7 k B T (middle) and -8k B T (bot- 
tom). Contours are drawn with 1 k B T spacing. 



above, the strands are modeled as 20 nm rigid rods. The 
only hybridization allowed is between strands on oppo- 
site plates, and we consider solution hybridization en- 
ergies of — 8/cb? 1 , —7ksT and — 6 k^T. Figure [5] shows 
the potentials of mean force of these plates as a func- 
tion of plate separation and lateral displacement (zero 
lateral displacement means that the centers of the circu- 
lar patches line up vertically). Such detailed potentials 
could be used to inform parameter choices for coarse- 
grained models of patchy particles^ when the patches 
are made from patterned DNA coatings. For nanoscopic 
colloids, we suggest that such patches could be realized 
by building the particles through DNA origami 35 . How- 
ever, to our knowledge, such precise experimental con- 
trol of coating patterns for microscopic colloids remains 
elusive, although recent experiments with triblock Janus 
particles are very encouraging^. 



B. Colloids coated with long ssDNA tethers 

To describe more flexible ssDNA tethers, we use 
the representation recently employed by Rogers and 
Crocker 17 to describe experimental measurements of 
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FIG. 6. Calculated number of bonds formed between the A- 
and B-type colloids described in Ref. 1171 when separated by 
5nm (top curve) to 30 nm (bottom curve) in steps of 5nm, 
averaged over 14 realizations of tether grafting points on the 
colloids. Inset: Same before averaging, i.e., for an individual 
realization of tether grafting points. 



DNACC interaction potentials. They model ssDNA teth- 
ers of contour length 40 nm as 8-segment freely jointed 
chains of Kuhn length 5 nm. The sticky ends that cap 
them are modeled as point-like, so that when two tethers 
bind, t hey form a 16-segment freely jointed chain. Ap- 



A2 



pendix 
such tether's 



describes how to calculate AG 



(cnf) 



for two 



i and j. 

In Rogers and Crocker's experiment, 
complementary types of sticky ends, A 



there are two 
and B, whose 

sequences are predicted^ by DINAmeltpSl to have bind- 
ing enthalpy AH° = — 47.8 kcal/mol and binding en- 
tropy AS = -139.6 cal/(mol-K). At temperature T, 
the A-B binding free energy in solution is taken to be 
AG = AH°-TAS°. Tethers with these sticky ends are 
coated on three types of spherical colloids, all of diame- 
ter 1.1 /im, A- type colloids are coated with 4800 ±480 A- 
type tethers, S-type colloids are coated with 4200 ± 420 
-B-type tethers, and AB-type colloids are coated with 
2400 ± 240 ,4-type and 3500 ± 350 B-type tethers. In all 
the following simulations, we have used the central value 
of all ranges. 

We first address the interaction between A and B col- 
loids. The coating densities used by Rogers and Crocker 
are sufficiently low that local grafting density fluctua- 
tions are relevant. We have thus prepared 14 indepen- 
dent configurations of uniformly coated A and B colloids. 
Figure [6] shows the number of bonds formed between two 
such colloids when separated by 5 nm to 30 nm versus 
sticky end binding strength AG , as calculated by ex- 
plicit Monte Carlo simulation and the explicit SCT the- 
ory, averaged over the 14 independent coating realiza- 
tions. The agreement is nearly quantitative, both when 
averaged over grafting realizations and independently for 
each such realization (inset). 

By using Eq. (15), we can integrate the curves in Fig- 
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FIG. 7. Potential of mean force between A and B colloids of 
Ref.[T7|at temperatures of 36.0°C, 33.0°C and 30.5°C (upper, 
middle and lower curves). See text for details of the various 
methods. 
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FIG. 8. Potential of mean force between two AB colloids of 
Ref.Q7|at temperatures of 32.0°C, 27.0°C and 24.0°C (upper, 
middle and lower curves). See text for details. 



urc|6]to obtain the mean attraction between the A and B 
colloids. Further accounting for the entropic repulsion 
caused by confining the tethers between the two col- 
loids results in a potential of mean force for this sys- 
tem, shown as dots in Figure [7] In this Figure, we also 
show the result of calculating the colloid interaction po- 
tential by calculating the interaction energy per unit area 
between two similarly coated plates, and then applying 
the Derjaguin approximation. The interaction free en- 
ergy between plates is calculated in three different ways: 
by explicit MC simulations, by using SCT, and by us- 
ing its mean-field approximation, SCTMF. Our results 
here help establish the Derjaguin approximation as valid 
for ssDNA tethers modeled as freely-jointed chains. For 
completeness, we show similar results for the AB-AB sys- 
tem explored by Rogers and Crocker in Figure [8] 

It should be noted that the simple freely-jointed chain 



model of ssDNA used here and in Ref. [T7] may not ad- 
equate ly ca pture the complications of the experimental 
systerrPSEll Our aim in this paper is limited to estab- 
lishing the SCT as an accurate and fast alternative to 
detailed MC calculations, regardless of the underlying 
model used to calculate the configurational entropic fac- 
tors associated with individual tethers binding. 



IV. SUMMARY OF RESULTS 

For ease of reference, in this section we summarize 
the main results of the paper. These results are im- 
plemented in the Python package DNACC, available at 
http : //github . com/patvarilly/DNACC. 



A. Self-consistent theory (SCT) 

Our central result describes the interaction details of a 
DNACC system where each strand i is modeled explicitly. 
The algorithm is: 



1 



Using a suitable microscopic model, calculate the 
binding free energy AGjj between every pair of 
strands i and j. When i equals j, or when i and j 
cannot bind, set AG^ to infinity. This binding free 
energy is the sum of two contributions: the solution 
hybridization free energy of the sticky ends, AG".,-, 
and the entropic cost of binding of the two teth- 
ers that hold the sticky ends, AG l -™ f ' 1 . We provide 
expressions for AG|™^ for rod-like and polymeric 



tethers in Appendix^A"} 



2. Solve the self-consistent set of equations for each 
of the quantities pi, i.e., the probability pi that 
strand i is unbound: 



1 + E, e-z^y 



for all 



Pj 



(33) 



It is straightforward to solve these equations using 
fixed-point iteration, i.e., making an initial guess 
for {pi}, then repeatedly inserting the values of 
{pi} into the right-hand side of the above equa- 
tions to obtain an improved estimate for {pi}. A 
solution can also be obtained through minimization 
of a function with no other local minima or saddle 
points, as detailed in Section [HA] 



Calculate py, the probability 
i and j are bound: 



Prj 



PiPj- 



that strands 



(34) 



The interaction free energy of a DNACC system can 
then be calculated through thermodynamic integration 
of the binding free energies AG^, as follows: 
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1. Calculate the the free energy of repulsion, F rep , 
equal to the free energy of the system when none 
of the strands can bind, i.e., when all the AG^- are 
infinite. This is simply the sum of the free energies 
of repulsion of the individual strands. We assume 
this quantity can be calculated analytically or oth- 
erwise, as described in Appendix [A"| 

2. Define a thermodynamic integration variable, A, 
and denote by Pij(X) the probability that strands 
i and j are bound when the all the binding free 
energies are respectively set to AGV, + AfceT, cal- 
culated as above. 

3. The free energy of the system is given by the fol- 
lowing equation 

/>oo 

F = F rcp -k B T / dAVV-(A). (35) 



B. Self-consistent mean-field theory (SCMFT) 



The free energy per unit area between two DNA-coated 
parallel plates separated by a distance h is estimated as 
follows: 



1 



Calculate the free energy density of repulsion, 
ftep{h), equal to the free energies of repulsion of the 
individual st rand types, weighted by their coverage 
density (Eq. ( |A8| ). Expressions for this quantity in 
the case of rod-like and polymeric ligands are given 
in Eqs. dA8b, (|A15|) and (1X331. 



2. Define a thermodynamic integration variable, A, 
and denote by <r a b(A) the number of a-b bonds per 
unit area when the all the solution binding free en- 
ergies are respectively set to AG° 6 + Ak B T, calcu- 
lated as above. 

3. The free energy density of the system is given by 
the following equation 



/>oo 

f(h)=f Iep (h)-k B T / dA 
Jo 



a<b 



CTab(A). 



(39) 



When colloids are coated uniformly with a small num- 
ber of different DNA species at high density, it is useful 
to model the colloid interaction in mean field. Assuming 
that the colloidal surfaces have radii of curvature signifi- 
cantly larger than the range of interaction of two strands, 
we may first calculate the interaction free energy of two 
parallel plates as a function of plate separation, h, then 
apply the Derjaguin approximation. 

The details of the interaction between two parallel 
DNA-coated plates separated by a distance h are esti- 
mated as follows: 

1. For each pair of DNA species, a and b, calculate a 
binding factor K a b as follows: 



d 2 r b exp[-/3AG { ; 



(enf) 
b 



(36) 



where AG^ 11 ^ is the configurational cost for a fixed 
a- type ligand to bind a b- type tether grafted at r&. 
If a DNA species is present in both parallel plates, 
treat the DNA on each plate as a separate species. 
Simple analytic expressions for K a b for rod-like and 



polymeric tethers are given in Eqs. (251 and (26) 



2. Solve the self-consistent set of equations for each 
of the quantities p a , i.e., the probability that an 
a-type strand is unbound: 

Pa = rrp4 , (37) 

1 + l^b K abVbPb 

where a\, is the number of 6-type strands per unit 
area. 

3. Calculate the number of a- and 6-type bonds per 
unit are, a a b, as follows: 



Finally, to estimate the free energy of interaction be- 
tween arbitrary colloids, use the Derjaguin approxima- 
tion^. I n the case of two spheres of radius R\ and i?2 , 
respectively, with a distance of closest approach h, this 
approximation yields: 



F(h) 



2tt 



R t 1 + R 2 1 Jh 



dh'f(h'). 



(40) 



Cab = &aPaKabPb&b- 



(38) 



In using the Derjaguin approximation, it is important to 
shift all free energies appropriately to ensure that f(h) is 
zero at infinite h. 



V. CONCLUSIONS AND OUTLOOK 

The theory for DNACC interactions that we have pre- 
sented, both in its explicit form (Eqs. |6|, f ty and (15)) 
and its mean-field approximation (Eqs. (|2l|)-([24|) is ca- 
pable of describing many different DNACC setups, vary- 
ing in grafting density, sticky end interactions and tether 
details. The theory is robust and correctly incorpo- 
rates the statistical mechanics of valence-limited binding, 
which plays a key role in the behavior of DNACCs, and 
it does so at much lower cost than full MC simulations. 
Extending the theory to treat general ligands within our 
approximation is straightforward, as exemplified in Ap- 
pendix [X] for freely-jointed-chains. 

When comparing our full theory to MC data in the 
context of DNACCs, we have found quantitative agree- 
ment. The mean-field approximation, which substan- 
tially generalizes earlier approaches, is in quantitative to 
semi-quantitative agreement with the exact results. Sig- 
nificant discrepancies arise only when the grafting density 
fluctuations that are neglected in the mean-field approach 
become relevant. 
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Our theory accounts for all the complications arising 
from valence limitation and combinatorial entropy. Nev- 
ertheless, the binding free energy AGy of a pair of ligands 
i and j is taken as an external input, and the agreement 
of the theory's predictions with experimental results will 
depend significantly on the quality of the estimates for 
the values of AG^. Simple models for generating these 
estimates, such as those used in this paper, can inform 
experiments, though apparently, they are not yet suffi- 
ciently accurate for quantitative agreement with experi- 
ments, as we have discussed in Ref. [T5J 



We expect our theory to be a significant aid in the 
search for viable approaches to design crystal structures 
a priori. Presently, existing free energy methods can 
establish the stable crystal structure of a system given 
its interaction potential. Conversely, robust methods 
are being developed for solving the inverse problem, de- 
termining the interaction potentials that might yield a 
particular target structured Designing colloids that re- 
alize such interaction potentials experimentally is an- 
other problem. In this regard, DNACCs, with their 
highly tunable interactions and vast parameter space 
(e.g., particle size, ligand structure and flexibility, coat- 
ing density, nucleotide sequence length and hybridiza- 
tion strength) , are a promising tool to build colloids with 
programmable interactions. For this reason, our theory 
can complement existing tools in crystal structure de- 
sign by quickly suggesting experimentally realizable con- 
structs that result in particular colloid-colloid interac- 
tion potentials. Whereas for simple potentials one can 
envision a manual search of parameter space guided by 
intuition and experience, it might also be possible to re- 
cast the problem in terms of a global minimization al- 
gorithm that optimizes the DNACC structures to obtain 
interaction potentials that are close to some target. It 
is precisely to help others easily use our tools in these 
and other ways that we have written a comprehensive 
and well documented series of Python modules that im- 
plement our theory, available under the GNU license at 
http : //github . com/pat var illy/DNACC . 



While we have developed our theory and the accom- 
panying programs with an eye towards applying them to 
DNA-coated colloids, other interesting applications are 
possible. In particular, our theory can help model any 
system where attraction is dominated by linker-receptor- 
type bonding. These include patchy particles, where the 
patches are realized by a nonuniform ligand coat. More 
generally, many biological systems could be treated with 
our theory, for example, interactions between cells (or be- 
tween cells and a substrate) that are mediated by mem- 
brane proteins. For these more general cases, care should 
be taken to include the effects of surface elasticity and 
how the entropic costs of ligand binding are, in turn, 
modified. 
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Appendix A: Entropic cost of ligand binding 

In this appendix we describe how to calculate the free- 
energy change AG^- of binding two ligands, i and j, 
when the binding interaction is physically limited to a 
small, reactive sticky end. We begin with a general 
treatment, and later specialize to the two important 
cases of short, stiff polymers and of freely-jointed chains. 
These two models have been used to describe the DNA 
strands on DNA-coate d colloids when these are short, 
rigid dsDNA E^ 15 1 16 1 25 1 and when they are long, flexible 
ssDNAl^li^. These models may also be useful to describe 
ligand-mediated binding in important biological systems, 
such as cells, where adhesion is, in many cases, medi- 
ated by different glycoproteins coated on their surface 
membrane^. 

Let Qi, Qj and Qij be the partition functions of the 
two unbound ligands and the bound complex that they 
form. In terms of these quantities, AGij is given by 

exp(-/3AG iJ ) = ^-. (Al) 

As in the main text, we have assumed that the unbound 
ligands do not otherwise interact with each other. The 
Boltzmann factor exp(— (3 AGij) is generally smaller than 
that of the sticky ends in solution at standard conditions, 
exp(— /3AG? ), owing to the entropic penalty of con- 
straining the configurations available to ligands i and j 
when their reactive sticky ends are bound. This defines 
the configurational free energy cost of binding, AG^ nf ^ , 
as follows: 

AGg nf) = AG iS - AG%. (A2) 

A realistic model that reproduces AGy for ligands 
such as DNA strands is expensive, as one would need 
to at least model each nucleotide of each strand explic- 
itljE^. It is computationally more convenient to separate 
the binding of the sticky ends from the configurational 
freedom of the rest of the ligand, and use a coarse-grained 
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representation of the sticky ends where AG^- is intro- 
duced as an explicit external parameter. Concretely, the 
internal partition functions of the free and bound sticky 
ends in solution, Zi, Zj and Zij, arc related to the equi- 
librium constant of binding^, 

exp(-/3AG° ) _ 2^ 



K 



Po 



(A3) 



where po is the standard concentration. Using this rela- 
tion in Eq. ( |Al[ ) to factorize the unknown dependence of 
Qij, Qi and Qj on the internal partition functions, we 
obtain 



Qi 



)(cnf) 
'ij 



QiQj 



(A4) 



l(cnf) 



)( cnf ) anc j Q^ ni ^ are approximately inde- 



where , ^ _ ^ 

pendent of the chemical identity of the sticky end ( "enf " 
stands for "configurational" ) . It follows that 



exp[-/?AG^] 



1 



j(cnf) 
'ij 



)( cnf )Q(cnf) ' 



PO i 



(A5) 



We stress that this factorization is an approximation that 
holds only when the binding of the sticky ends is only 
weakly affected by the presence of the rest of the ligand, 
as occurs when the ligand is an ideal polymer or when 
the binding is very strong. However, it is a useful ap- 
proximation because once AGg nf) has been measured or 
calculated, one can predict the values of AG^- for families 
of ligands that differ only in their sticky ends, providing 
us with a powerful design tool. 

We now specialize to the case of polymeric ligands with 
one end fixed and the other end sticky, such as polymers 
grafted at a surface. The main complication in treat- 
ing polymers is dealing with excluded volumes. For in- 
stance, in DNACCs, the DNA strands cannot penetrate 
the colloids. In the discussion below, we denote by r< the 
grafting point of polymer i. Given two polymers i and j, 
we denote by the vector r, — with length r^. We 
denote by £1 the region of space that can be explored by 
the monomers of the polymers. Polymers restricted to 
f2 are deemed "confined", whereas those for which this 
restriction is absent are "unconfined" . 

We note that the free energy of repulsion due to a poly- 
mer i, denoted by F rePj i, is captured by the dependence 

of Q^ cnf ^ on the confinement region Q. Indeed, letting 
Q( cnf ) k e p ar tition function of polymer i in the ab- 
sence of confinement, we have 



g-Wop,, =Q (cnf) / g(cnf)_ 



(A6) 



From this expression follows the total free energy of re- 
pulsion, 



,F le 



(A7) 



Similarly, the mean-field free energy density of repulsion 
is given by 



(A8) 



where a a is the coverage density of species a and F Iep>a is 
the free energy of repulsion of any representative et-type 
strand. 



1. Short, stiff polymers as ligands 

A short, stiff polymer, such as a short piece of dsDNA, 
can be modeled as a rigid rod, as has been done previ- 
ously in Refs. |T5J [TH] and [S5J We gather the results here 
for completeness and to clarify the origin of a Jacobian 
prefactor absent in Ref . fl5l but addressed in Ref. [TBI For 
rigid rods, the configurational partition functions are: 



Qg nf) = 


/ dr<5(|r - r 4 | 
Jn 


-Li)5(\ 


Q (cnf) = 


/ dr5( r-ri| 
Jn 


-Li), 


g (cnf) = 


/ dr<S(|r-r 3 | 


-Lj), 



(A9b) 
(A9c) 



where Li and Lj are the lengths of rods i and j. 

The first integral is most easily evaluated in cylindrical 
coordinates, with the z axis parallel to r^. Special care 
must be taken to isolate the Jacobian factor that arises 
when dealing with the product of the two (5-functions. 
After evaluating the integrals, we obtain, 

Q( cnf ) _ i - 



sin( 7l +7j)' 
Qf ni) = Qj. 



(A10) 

(All) 
(A12) 



Here, 0^- is the portion within Q of the circumference 
of the circle explored by the common endpoint of the 
hybridized rods, Qi and fij are the areas within f2 of the 
sphere that the endpoints of the free rods can explore, 
and 7^ and "/j are the angles that the hybridized rods 
make with ry. The factor of 1/ sin(7j +7v) was obtained 
in Ref. QI)] using two different representations of the 5- 
function, but the procedure here highlights the generality 
of the result. 

When the colloids are plates, one at z = and one at 
z = h, the region 17 is simply that where < z < h, 
and one can write an explicit expression for fi^-. For 
the case where the rods are grafted on opposite plates, 
such an expression is given in Appendix A of Ref. [TH] (see 
also Figure 8 therein). When the rods are grafted to the 
same plate, f^- is zero unless \Li — Lj\ < rij < Li + Lj, 
in which case 



LiSin^i I 2 arcsin[/i/(Li sin7i 



h > Li sin7j, 
h < Li sin 7j . 

(A13) 
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The quantity 0, also has a compact expression, 
r 2nLf, h>L h 



27rii/i, < h < Li. 



(A14) 



The free energy of repulsion of a strand i follows imme- 
diately from Eq. ( A6 ) , 



rep, i 



0. 



-knTln(h/Li), < h < U 
h > Li. 



2. Long, flexible polymers as ligands 

Following Refs.[17]and[T8l we can model a polymer i as 
a freely-jointed chain (FJC) with segments of length £. 
Let rj jS be the endpoint of segment s, with 1 < s < N. 
For convenience, let r^o = r i, and let Ui s = rj g — Ti, a —i- 
The chain segments are confined to the region f2, but do 
not otherwise interact. It follows that 



l(cnf) 



J|du M 5(|u 4 



£). 



(A16) 



where, abusing notation somewhat, the integration do- 
main confines the polymer to the region tt. A completely 
analogous treatment yields Q^"^ for a FJC with Nj seg- 
ments of length I, grafted at Yj . 

One way of modeling the binding of i and j is to co- 
alesce the endpoints of the two chains, so that r^jvi = 
Yj : Mj , resulting in a FJC with endpoints fixed at and Tj 
and with iVhyb = N$ + Nj segments {u s }, where 1 < s < 
iVhyb- However, if the persistence length of the polymers 
is short, then the binding site may have a size comparable 
to the segment length i. This is the case for ssDNA poly- 
mers with sticky ends that are 6 to 12 nucleotides long. 
Hence, we also consider a more general model for binding 
of i and j in which the final segment of the two polymers 
coalesces, so that r^jy, = r^-jy _i and jy 4 _x = Tj,Nj, 
resulting in a FJC of N hyh = N t + Nj - 1 segments^. 
The partition function of the bound system is 



(A17) 



„ Nhyb ^hyb 

^S nf) = / [II du s <5(|u s | -ejlsfa - £ u s ) , 



(A18) 



where is the partition function of a confined FJC 



' >3 



with fixed end points, and r^. The additional factor 

in Q-™^ results from fusing the end segments of the in- 
dividual polymers. 



Although the integrals in Eqs. (A16| and (A18) cannot 



be evaluated analytically even in simple geometries, as 
was done in the previous section, one can calculate their 
value numerically using standard polymer techniques. 



Let Qi(N) = (4tt£) n be the partition function of an un- 
confined ./V-segment FJC fixed at one endpoint, and let 
Qi (N, R) be the partition function of the same FJC with 
fixed end-to-end distance R. Their ratio yields the prob- 
ability density for the end-to-end distance of the FJC to 
be R, 



p(N,R) 



i(N,R) 
UN) ' 



(A19) 



which is known in closed forrr l 44 * 45 4 Using these partition 



functions, we can write Qij nf ' in the following convenient 
form: 



Q 



(enf) 



(enf) 



Q2(N h y h ,rij). 



(A20) 



The ratio on the right-hand side of Eq. (A20) is the 



probability that a random walk from to Vj lies wholly 
within £1 (i.e., does not penetrate the colloid). To com- 
pute it, we use Rosenbluth sampling^, a common tech- 
nique to generate a properly weighted ensemble of such 
FJCs. Specifically, to grow each segment s of one chain 
that starts at r*j, we first generate K possible segments 
U-Sj (1 < k < K), sampled with probability density 
p(Nhyb — (s + 1), \vj — r s+ i|) using rejection sampling 46 . 
We then pick at random one of the m s segments in 
{uSj} wholly within fl. The Rosenbluth weight of that 
chain is defined as 



m s 
~K 



(A21) 



The fraction of FJCs that lie wholly within is then the 
average Rosenbluth factor, 



;(cnf) 



i(iVhyb,7y-) 



(W) 



(A22) 



where the subscript implies that the average is computed 
as described above. 

We can estimate Q^ cnf ^ and Qj^ analogously, except 
that at each step in growing the chain, now of length iVj 
or Nj, the possible segments {u^i} are chosen uniformly 
on the surface of a sphere of radius I, since the end of 
the chain is free. We obtain 



Oi(JV ( ) 
Qf ni) 
Qi(Nj) 



(W)i 



(W) 3 



(A23) 
(A24) 



Substituting Eqs. (A17|-(A24| into Eq. (A5|, we ob- 



tain the following expression for the configurational en- 
tropy cost of binding i and j: 



exp(-/3AG'^) = p(iVhyb ' ry ' ) {W) ^ 



Po (W)i(W)j 



(A25) 
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FIG. 9. AG'™ f ^ for a plane geometry and the polymeric 
constructs of Sec. |III B| as obtained using MC simulations 
(symbols) and a harmonic spring approximation (Eqs. ( |A27[ ), 
( A29| and JA30I), full lines). The two FJCs are tethered either 



to different plates (above) or to the same plate (below). 

In Figure [9] we show the configurational entropy cost 
of hybridizing two FJCs (made of 8 segments with £ = 
5 nm) tethered to parallel plates. Interestingly, if AG 
is expressed in terms of the distance between tethering 
points, Tij, and the plate separation, h, the dependence 
on the latter is marginal. This suggests that, for ideal 
chains, the cost of stretching a chain (p(7Vh y b, ry) in 
Eq. ( |A25l )) is more important than the excluded volume 

terms (W). Hence, we may approximate AGj™^ analyt- 
ically by neglecting all chain degrees of freedom but the 
sticky end. In particular, using a Gaussian approxima- 
tion for p(N, rf^, 



p(N, r) 



3 

2nN£ 2 



3/2 



exp 



-3r 2 
2NP 



(A26) 



we can model the ligand as a reactive point-like sticky end 
that is confined between the two walls, and is attached 
to the tethering point via a harmonic spring with spring 
constant k — "ik^T/Nt 2 . The partition function of an 



unhybridized polymer is then given by 



Q| cnf) « f Avp{N,\v-v l \) = \ C ii 
Jn 1 



3 

2NP 



, (A27) 



where the integration domain, SI, is [-co : oo] 2 x [0 : h], 
and the error function, erf(x), is defined as in Eq. ( |27[ ). 
Similarly, the partition function of the hybridized poly- 
mer is given 



Q4 cnf) « J dr p(N, |r - ri |) • p(N, |r - rj |). 



(A28) 



When i and j are tethered to opposite plates, so that the 
bound pair forms a bridge, we obtain 



Q 



(bridge) 
ij 



1 

8 V ttN£ 2 



3/2 



erf 



: 



3r 2 



h 3 

2\lN£ 2 ) eXP [ 4N£ 2 
(A29) 

whereas when i and j are tethered to the same plate, so 
that the bound pair forms a loop, we obtain 



Q 



(loop) 



1 



16 \ttN£ 2 



3/2 




rf h 



exp 



3r 2 



AN£ 2 
(A30) 

From these formulas, using Eq. ( A5 1 , we obtain analytical 

estimates of AG^"^. In Figure jo| these estimates are 
shown in the solid curve. We see an overall agreement 
with the exact result up to = 50 nm, where the finite 
extensibility of the chain becomes important. 

The above formulas can be used to derive an approxi- 
mate expression for the coefficients K a t used in the mean- 
field treatment of Section II C (Eq. (26 1). Indeed, for 
bridges, 



and for loops, 



^(bridge) 

dr a6 27rr 6^ c ^ (cnf) . 



dr ab 2irr ab 



q1T p) 

Oi cnf) Qi cnf) ' 



(A31) 



(A32) 



In both of these cases, we have assumed that Q^ b ^ goes 
to zero sufficiently quickly for the integrals to be insensi- 
tive to their upper bounds. 

Finally, the above formulas also yield simple approxi- 
mate expressions for the free energy of repulsion, -Frep,i, 
of a strand i. From Eq. (A6), we obtain 



rep, i 



-/cB^ln erf 



ih 2 



2N£ 2 



(A33) 



Appendix B: Relation to Local Chemical Equilibrium 
treatment 

The derivation of Eqs. ^ and ([7| can be regarded 
as a "ligand-centric" result: the focus is solely on which 
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ligands (DNA strands) bind to one another. We can take 
a complementary "sticky-end-centric" approach to obtain 
the same results, which sheds light on the recent Local 
Chemical Equilibrium treatment of Rogers and Crockei^ 
and corrects its main deficiency. 

As in the main text, we first treat all strands sepa- 
rately, then make continuum approximations for simpli- 
fied geometries. Following the notation of Ref. [T71 let 
G°(r) be the normalized probability density of finding the 
sticky end of tether i at r, before any hybridization oc- 
curs. For this treatment only, we assume that the strands 
interact "ideally", i.e., with an interaction energy AG°j 
when their sticky ends coincide and none otherwise, so 
that the configurations of the tether of i with sticky end 
at r have the same statistical weight whether the tethers 
are bound or unbound. In terms of partition functions, 
strands i and j interact ideally if their partition function 
when bound, Qij, is given by 



Qii = / drQ i (r)Q,(r) 



ex P (-/3AG°.) 



Po 



(Bl) 



where Qi (r) is the partition function of unbound strand i 
with sticky end at r, and Qj(r) is defined similarly. The 
restriction of ideal binding, absent in the general treat- 
ment of the main text, is satisfied when the tethers are 
ideal chains, for example. With this additional assump- 
tion, we rewrite Eq. ([I]) as follows: 



Z 



e n 



drG°(r)G°(r) 



cxp(-/3AG 



Po 



(B2) 



where po is the standard concentration, 1M. The analog 
of Eq. ^ is: 



We can then finally write the self-consistent condition for 
the quantities {pi} as follows: 

1 = Jdr /^(rH^G^r)!, (B6 



G i (r)=p i G°(r). 



a) 
(B6b) 



Comparing the above condition to the equations in 
Ref. [TTJ we conclude that the "local chemical equilib- 
rium" approximation introduced by the authors is equiv- 
alent to our approximation Z-^-j/Z « PiPj, but that the 
self-consistent condition is significantly different from the 
condition that they proposed, 



G(r) = G°(r)- YMt). 



(B7) 



Physically, in the treatment of Rogers and Crocker, de- 
pletion of free sticky ends i at a point r establishes a den- 
sity gradient in Gj(r) beyond that implied by the poly- 
mer statistics of tether i, but this density gradient is not 
allowed to relax, effectively imposing a nonequilibrium 
condition. Our treatment, summarized in Eq. (B6), cor- 
rects this shortcoming. 

As in the main text, for those cases where there is 
translational invariance, we may approximately treat to- 
gether tethers with similar binding behavior. In the con- 
crete case of two plates perpendicular to the z axis, and 
in the absence of binding, the concentration of a-type 
sticky ends at a height z is: 



G a °(z) = a a P a (z), 



(B8) 



where 



= / dr ^C°(r) 



Ec?(r)c .( r) f5h^%i> , mi 



Po 



As in the main text, let pi = Z-i/Z, and approximate 
Z-i-j/Z by PiPj. As in Ref. [T71 let Cj(r) be the con- 
centration of free sticky end i at r, equal to piC°(r) by 
comparison to the equation above. There being only one 
such sticky end, Ci (r) is also the probability density that 
the sticky end of tether i is free and is at r. We obtain 



1= / dr C.W+^QWQfr 



exp(-/3AG?- 
Po 

'(B4) 

More concisely, let Cij (r) be the concentration of bound 
sticky ends i and j, given by 

exp(-/3AG°) 

C,s(r) = ^(r C,.(r)-^: *lL. (B5) 

Po 



P a (z)= / dxdyCf(x,y,z), 



(B9) 



and i is any representative tether of type a. The corre- 
sponding mean field self-consistent equations are 

°a =Jdz ^C a (z)+J2 C ab(z)^ , (BIO) 



C a (z) = Pa C° a (z), 



(Bll) 



G Q b {z) = C a (z)C b (z) eM /?AG ° b) . (B12) 

Po 

Rearranging these equations to obtain explicit equations 
for {p a } yields 

Pa = 7~i r j n , in / ,cx P (-/iAG , ) ( B13 ) 

1+PbPb J dzP a (z)P b (z) PV ' pQ ab! 



Comparing this last equation to Eq. (21 1, we find that 

-|0 



-/3AG u a 



PO 



dzP a (z)P b (z). 



(B14) 
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This is a specialization of Eq. ( 22 ) under the additional 



assumption of ideal binding interactions. 

For the concrete case of short, rigid- rod-like tethers of 
length L, the probabilities P a (z) are particularly simple, 
and depend only on rod length, plate separation, and the 
plate on which a-type tethers are grafted. In particular, 
if they are grafted on the plate at z = and the other 
plate is at z — h, then 



Pa(z) = { 



l/h, < h< L] 

l/L, h > L and < z < L; (B15) 
0, otherwise. 



Otherwise, if a-type tethers are grafted on the z = h 
plate, then 



Pa{z) 




0<h< L; 
h > L and h 
otherwise. 



L < z < L: 



(B16) 



From these results, one can easily obtain Eq. (25 1. This 



derivation should be compared with the explicit integrals 
involved in calculating K a \> using Eq. (22) and the rele- 



vant expressions for AG^™^ given in AppendixlAlof this 
paper and Appendix A of Ref. HU We have verified that 
both routes give identical results when evaluated numer- 
ically for Lb — L a and Lb = 1.7 L a . 

In the absence of translational invariance or some sim- 
ilar symmetry, it is nontrivial, in general, to write down 
self-consistency conditions similar to Eq. ( B6 ) where like- 



type tethers are grouped together. Unlike the incorrect 
condition of Eq. (B7), the correct self-consistency con- 



ditions are inherently non-local from the point of view 
of sticky ends: the change upon binding in sticky end 
concentration of a-type sticky ends at r depends on the 
precise identity of the a-type tether to which that sticky 
end is attached. Indee d, o ne needs to retain the level of 
detail inherent in Eq. (B6). As such, there is no advan- 
tage in adopting a "sticky-end-centric" point of view as 
opposed to a "tether-centric" point of view when deal- 
ing with arbitrarily shaped colloids, polymer statistics 
or non-uniform coverage. Moreover, the tether-centric 
equations are simpler, since they result from integrating 
over the irrelevant degrees of freedom describing the po- 
sition of the sticky ends, and are more general, since they 
do not assume that binding is ideal. 

To conclude this section, we present a numerical com- 
parison between Monte Carlo results, our self-consistent 
theory and the LCE theory that illustrates the extent to 
which these approximate theories deviate from an exact 
treatment. We consider two periodic parallel plates, each 
with 500 tethers, at a density comparable to that used in 
Ref. E] (4500/ (4tt x 550 2 ) = 0.001183 tethers/nm 2 ). The 
tethers are modeled as ideal 8-segment chains of Kuhn 
length 5 nm that, when bound, form a 16-segment ideal 
chain, as described in Appendix [A] Figure [lO] shows the 
number of bonds formed in this system as a function of 
the solution hybridization free energy, AG , of the sticky 



Self-Consistent Theory 
Local Chemical Equilibrium 
Monte Carlo 




/3AG 



FIG. 10. Exact and estimated number of bonds formed be- 
tween two ssDNA-coated parallel plates separated by 15 nm 
(lower curves) and 25 nm (upper curves). The Self-Consistent 
Theory is that of Eqs. (J©J) and Q, while the Local Chemical 
Equilibrium curve results from the methods of Ref. 1171 See 
text for details. 



ends of tethers on opposite plates, when the plates are 
separated by 15 nm and by 25 nm. As is observed, our 
self-consistent theory closely tracks the Monte Carlo re- 
sults, whereas the estimates obtained using the LCE ap- 
proach of Rogers and Crocker differ appreciably from the 
exact ones. 



Appendix C: Self-consistent equations from a Saddle Point 
Approximation 

In this appendix, we present an alternate derivation for 
the self-consistent mean field theory in Eqs. ( 21 ) and ( 22 ) 



in two specific setups, based on a saddle point approxima- 
tion of a mean-field partition function. The two geome- 
tries we explore are the first two examples in Section [III| 
which we refer to as the single linkage system and the 
competing linkages system. 



1. Single linkage system 

In the single linkage case, two parallel plates of area A 
are respectively decorated with N a and Nb strands of 
two complementary kinds of DNA, a and 6, whose sur- 
face concentration is thus a a = N a /A and ab = Nb/A. 
Without loss of generality, we assume that N a < Nb- Let 
n b<-a be the average number of o-type strands that a 
fixed a-type strand can bind to. When no bonds have 
yet formed rib^ a = VbAb^a, where Ab^ a is the area 
that encloses all possible grafting points of &-type bind- 
ing partners of a fixed a-type strand. However, after 
having formed n a b bonds between a- and 6-type strands 
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the average number of available linkers is reduced tc5^ 



N b - n ab 

nb^a{n a b) = 77 VbA b< 

Mb 



(CI) 



The average Boltzmann factor for binding a- and 6-type 
strands is K a b/Ab^ a , where K ab is defined in Eq. (22 1. 



With this results and the previous one, we can write a 
mean-held partition function of the system as follows 



E 



n a b=0 



N a 
n ab 



Ha!,— 1 js 

II a(i) 



L j=0 



A,, 



(02) 



We now isolate the dependence of Z on A in order to 
take the limit A — > 00, which we then evaluate through 
a saddle-point approximation. Let s ab — n a b/A, and 
approximate the sum in Eq. ( C2 ) by an integral with 



respect to s ab . Assuming the maximum of the integrand 
is far from the integration limits, we obtain 



Z« A 



J ds ab 



exp[A -S(s ab )], 



(C3) 



where, using Stirling's approximation, 

S(s) — cr a In a a + (Tb In a b — s In s — s 

- (a a - s) ln(a a - s) - (a b - s) ln(a b - s) 

+ slnK ab . (C4) 



In the limit A — > 00, the partition function will be dom- 
inated by the values of s that maximize S(s). Setting 



dS/ds to at s = a ab , we obtain the following equation 
for the average density of hybridized strands, a a y. 



(Tab = K — Vab)K a b{<Jb ~ CT a b), 



(C5) 



which restates Eq. ( 23 ) after identifying a a p a with a a 



o a b and likewise for <JbPb- Using Eq. (21 1 one can also 



rewrite Eq. ( C5 ) as 



Oab 



a a K a b{ab — &ab) 



(C6) 



1 + K ab (a b - a ab ) ' 
which is the expression used in Ref. 1161 

2. Competing linkages system 

In this section we generalize the saddle point approach 
to the competing linkages system. Like in Section III 



we have two parallel plates of area A, one grafted with a- 
and 6-type strands, the other with a'- and fo'-type strands. 
The only allowed interactions are a-a' , a—b and a'—b' . An 
equivalent treatment applies when, as in Ref. 1161 allowed 
interactions are a—a', a-b' and a'—b, so that loops (intra- 
particle bonds) are forbidden. The only differences be- 
tween the two systems are the labels in the equations and 
the values of the K a b factors below. 

Using the same notation as in the previous section, we 
hnd that the average number of binding partners avail- 
able to a hxed strand of each type is 



T^a' <— a (j^aa' ? ^a'b' ) 



Ng> - Ugg, - Jig 

N a , 



N b - n ab 

nb^a{n a b) = 77 CT b Ab^ a , 

iVf, 



&a' A a f ^_ a , 

(C7a) 
(C7b) 



Hy «- a' [n a ,y ) = — ay Ay <_ „/ . (C7c) 

i\y 

Assume for simplicity that N a < N a >, N a < Nb and 
N a > < Ny ■ A nearly identical argument applies for a 
different relative populations of strands. A mean-held 
partition function for the system is 



N a N a—n aa ' N a ,—n aa / 

z= E E E 

n nn i— n ab — n n / h /—0 



N a 



n 

k=0 



K aa , 

A a i <_ a 



n a , ^ a (k,0) 



n ab 

N a N a -n aa , N a t—n aa i 

E E E 

n/ =0 n a h=0 n n / h /—0 



n ab — l 

n 

, t=0 



.4 



K a b ,.s 

n b ^ a {t) 



n a >y 



II A a ' b ' ny^a'iv) 



A 



riaa'J \n a a'J N, aa 



N a - n aa '\ ( N b \ n ab \ , 

— ^ {KabCTb T ab 

riab J \n a bj N h 



V "' " """'^ A ' V > " a ' V h Kal yOy) n ^ ■ (C8) 



n a >y J\na>yjN, 
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As in the previous section, we isolate the dependence 
of Z on A to take the A — > oo limit. Let s xy — n xy /A, 
and approximate the sums in Eq. ( C8 ) by integrals with 



respect to s aa ', s ab and s a ' b '- Assuming the maximum of 
the integrand is far from the integration limits, we obtain 

Z ~ A 3 J J J ds aa > ds ab ds a >b> exp[A ■ S(s aa , , Sab, Sa'b')}, 

(C9) 

where 

S(s aa ',S ab ,S a > V ) = COnst. 

+ s aa ' (In — — - 1) + s a6 (ln - 1) + s a > b > (In — — - 1) 

Saa' $ab $a'b' 

- {<Jb - s ab ) ln[cr b - s ab ] - {ay - s a > b >) ln[av - s a > b >] 

~ (<r a - $aa> - S ab ) ln[(7 a - S aa > ~ S ab ] 

- {<r a > - s aa > - s a > b >) ln[er a / - s ao > - s a / fe /] , (CIO) 

where the first, constant term does not depend on s xy . 

The point (a aa ' > 0q6, 0V&') at which S is maximized, 
found by setting all its partial derivatives to zero, yields 
the average bond densities for each pair of species. They 
satisfy 

Oaa' = (&a — a aa' ~ ^ab)K a a'{^a' — V aa' — Va'b') 

= <T a PaKaa'Pa'<Ja', (Clla) 
0"afc = (0"a - °aa' - ^ab')K ab {<T b - O a b) 

= (T a p a K ab p b o- b , (Cllb) 

O'a'b' = ("V - <Jaa' ~ Oa'b' )K a 'b' [pv ~ Oa'b') 

= <T a 'Pa'K a ,b'Pb"7b> ■ (Cllc) 

These equations correspond to the mean-field self- 
consistent theory of Eqs. (21) and (23)) for the com- 
peting linkages model. 

We now show that the saddle points equations ( |C11[ ) 
are also equivalent to the expressions derived in Ref. 1161 
(once a a b and o a 'b' ar e repla ced with a a b> and a a 'b re- 
spectively). Using Eq. (Cllc) two times, we have 



a a' — Vaa' — a a'b' 

OV — (T aa > 



a a 'b' 



(oV - a aa > - Oa'b') + &a'b' K a i b i[a b > — (Ja'b 1 ) 

1 + K a i b i{a b i — <r a ' b i) 



(C12) 



Substituting this result into Eq. (|Clla|) , together with 
the definitio; 

®aa' &a-K-aa' 



the definition of p a in Eq. (21 1, wc find that 

Cq' — <Jaa' 



1 + K a i b i{a b i — cr a > b ') 

1 _l_ Kaa'(Pa' ~ CT a a') \ V ( - „ \ 
1 + i^a'b'yCTb' — CT a 'b' ) 

&aK aa >(<7 a ' — &aa') ' 1 + K a ' b '{(Tb' - Pa'b') 

+ K aa '{a a ' — a aa ') + K a b(a b — a ab ) 

1 -l 

+ K ab (a b — <x ab ) ■ K a ' b i{a b i — a a i b t) 



(C13) 



This is exactly the expression used in Ref. US] to compute 
the density of a-a' linkages, derived there by applying an 
heuristic self-consistent prescription to an approximate 
estimate of a aa ' ■ A similar argument recovers the ex- 
pressions for a a b and o a 'b' used in Ref. [TBI 
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